Two-phonon states in alkali-metal clusters. 

Abdellatif Abada* and Dominique Vautherin^ 

Division de Physique Theorique, * 
Institut de Physique Nucleaire, 
F-91406 , Orsay Cedex, France 

Abstract 

Two phonon-states of alkali- metal clusters (treated as jellium spheres) are calculated 
by using a method based on a perturbative construction of periodic orbits of the time- 
dependent mean-field equations. Collective vibrations with various multipolar ities in 
charged Na^i are considered. 

PACS numbers : 36.40.+d, 31.50.+W, 33.20.kf 



Submitted to Phys. Rev. B 



IPNO/TH 92-80 December 1992 

* email : "ABADA@ipncls.in2p3.fr" 
Wail : "VAUTHERIN@ipncls.in2p3.fr" 

+Unite de Recherche des Universites Paris XI et Paris VI associee au C.N.R.S 



1 



2 



INTRODUCTION 

For metallic clusters, it appears reasonable to approximate the valence electrons as a system 
of independent fermions enclosed in a box, in close analogy with the mean-field description of 
nucleonic motion in nuclei. The many body methods developped for collective modes in nuclei 
are thus a natural tool to investigate, in a first stage, collective motion in metallic clusters. 
The length and the energy scales in metallic clusters and nuclei are of course widely different. 
The interactions also have a different nature. The strong force mediated by pion exchanges 
ensures nuclear binding while electrons interact with each other and with the nuclei through 
electromagnetic forces. 

The analogies between nuclei and metallic clusters became apparent in 1984 with the ex- 
perimental discovery of shell structure and magic numbers in sodium clusters by W. D. Knight 
et al [lj] and parallel theoretical efforts by W. Ekardt ||. Since then, the physics of metallic 
clusters has been of continuing interest ||. The valence electrons in these clusters can be de- 
scribed as delocalized particles moving in a mean-field due to their mutual interaction and to a 
uniform positive background. This approximation is known as the jellium model. As a result, 
the ionic background contributes only an electrostatic attraction to the mean-field, whereas the 
valence electrons are treated as quantum particles exhibiting exchange and correlation effects in 
addition to their electrostatic repulsion. This treatment is usually carried out in the framework 
of the density-functional formalism |5| . 

The experimental discovery of a giant dipole resonance in metallic clusters interpreted as 
a collective electronic excitation (surface plasmon) [7], |8| induced a substantial theoretical 
effort. To understand the presence of collective motion, the usual random-phase approximation 
approach (RPA) and the time-dependent local-density approximation (TDLDA) have been 
applied [9-14]. Our aim in this article is to describe collective motion of large amplitude 
in metallic clusters with a method which is more general than RPA. It is based on secular 
perturbation theory. More precisely we start from a perturbative method to construct finite- 
amplitude periodic orbits of the nonlinear time-dependent mean-field equations. This method 
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was developped recently to describe giant collective resonances in nuclear physics flTJ, [L6| and 
the Skyrmion breathing mode in hadronic physics fl7 |. Although this method is not always 



applicable due to the possible presence of branches and bifurcations |Tj| we found it well adapted 
to investigate giant collective oscillations in the case of nuclei. From the knowledge of periodic 
orbits we were able to calculate collective energy spectra by a semiclassical quantization or by 
a construction of the collective Bohr-type Hamiltonian. 

The present article is organized as follows. Sec. I presents a summary of the self-consistent 
spherical jellium-background model and defines our notation. In Sec. II we review briefly the 



expansion method of Refs. 16], [T^J and we specialize to the alkali-metal clusters case. In 
Sec. Ill we present our results concerning monopole, dipole, quadrupole and octupole modes 
in charged Nd^i (closed-shell spherical cluster with 20 electrons), and perform a comparison 
with available data JL9|, ^(J and other theoretical results ||, [11], [21]. In Sec. IV we discuss the 
applicability of the method to metallic clusters and summarize our main conclusions. 

1 The jellium model 

In the jellium model, the N ions are replaced by a uniform background of positive charge 
having radius R = r s N^ where r s = 3.93 a.u. for sodium. In the density functional approach 
the electrons gas energy E[p] is assumed to be a function of the one-body electron density 
matrix p(t) 

(r|p(t)|r') = E^(r,t)0*(r',t), (1.1) 

where the 0j's are the single particle wave functions of the A occupied single electrons orbitals. 
In this work, we use the density functional (see Ref. ||] for a general discussion ) 

E\p\ = /d 3 r { ir(r,*) - 2(I)V3 p V» (M) + P (r,t)£ xc (p) + \ / dV + 
Vj(r)p(r,t) } +Ejj . 

(1.2) 
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We use atomic units (e = m e — % — c = 1) for which energies are expressed in units of 
2R y = 27.2 eV (R y being the Rydberg constant) and lengths in units of the Bohr radius 
a = 0.53A. In Eq. (|1.2|) the first term is the kinetic energy density r(r,t) = YhLi |V0j(r, t)\ 2 , 
the second term is the Coulomb exchange in the local density approximation. The third one 
simulates correlations, while the fourth one is the direct Coulomb term. The fifth term is the 
Coulomb interaction of the electrons with the jellium. For an uniform density distribution for 
the N ions, the jellium coulomb energy Ejj is ^N 2 /R and Vj is given by 

f -N/r r > R 

Vj(r) = \ 3N ,,2 (1.3) 

For the correlation term, we follow the Gunnarsson-Lundqvist parametrization || : 

*~ = -°- 0333 G iik { i-/ I3] (L4) 

where the function G is defined by 

G(x) = (1 + x 3 )Ln[l + 1/x] -x 2 + x/2 - 1/3 . 



The time-dependent mean- field equation reads [22 



ip=[W(t),p(t)], (1.5) 
where W(t) is the mean-field Hamiltonian SE/Sp. It reads — |A + U(r,t) with 

U(r, t) = -(^) 1 /3 p i/3 (r . ; t) _ 0333 Ln[1 + n. 4 (^p(r, t)) 1 / 3 ] + / dV f^4 + V A r ) ( L6 ) 
7r 6 J |r — r I 

Finally, we recall that the one-body density matrix defined in Eq. (|1 . 1| ) must satisfy at all 

times the conditions Trp(t) = A (A being the number of electrons) and p 2 (t) = p(t) ||22|| . 

At this stage a word of caution about the application of the variational procedure should 

be given. The energy functional should be minimized with respect to the electron density 

only. In contrast it should not be minimized with respect to the jellium parameters. Indeed, 

such a procedure would lead to the prediction that potassium and sodium clusters have the 

same radius, which is in conflict with observations. This difference also suggests that it is not 

legitimate to perform a minimization with respect to the shape of the jellium. 
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2 The method 

The method we use to construct periodic orbits of the mean-field equation (|1.5|) with a given 
period T consists in expanding the one-body density matrix p(t) as a power series in the 
amplitude of the vibration e [TjJ. In order to take into account the dependence of the frequency 



uj = 2tt/T on the amplitude f23|, we expand p(t) as 



pit) = p + epi(— t) + e 2 p 2 (— t) + e 3 p 3 (— t) + ... (2.1) 

UJq LOq UJq 

with 

lu = ujq + tuJi + e 2 U2 + ... (2.2) 



The first order equation is solved by the usual ansatz p2 



pi(t) = r/exp (— iuot) + r/ + exp (+iu t) . (2.3) 



After inserting the expression of pit) given above in Eq. ( |1.5|) , we obtain the evolution equations 
of the successive terms pk (k = 0, 1, 2, ..) by identifying the various powers of e : 

[ Pk = {i W Ji Pm] -XljPm} ■ (2.4) 

j+m=k 

In Eq. (|2.4|) , the operator H^- corresponds to the j'-th order contribution to the mean-field 
Hamiltonian W(p) = J2k 6 while Qj is the j-th order contribution to the following expression 
1 _ too _ g _ ^ e l+l 'uoiuJv + ^ e w ' +r 'umup 



1=1 ^0 ^0 ^0 



The static mean-field solution po of Eq- fll.SQ corresponds to A; = in the notations of Eq. 
It satisfies the static equation [Wo, po] = where Wo is the static mean-field Hamiltonian. We 
define hole (h) and particle (p) states (normalized to one) by 

Po|<M = \<f>h), Po|0 P ) = 0, W \(f) h ) = e h \<f> h ), W \(j> p ) = e p \(f) p ) . 

In first order one recovers the well known linearized mean- field (RPA) equation ipi = M.p\ 
which is solved by the ansatz ([2.3|). The RPA matrix M. is defined by 

M Pl = [Wo, pr] + [Wi,p ] . 
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Since pit) must correspond to a Slater determinant, it must satisfy the condition p 2 (t) = pit) 
to each order in e. In fc-th order (k = 0,1,2,..), this condition reads 

Pk = PoPm • (2.5) 

j+m=k 

An immediate consequence of this last equation is that it provides at once the particle-particle 
and hole-hole matrix elements of the operator p^ from the knowledge of the lower order terms 
Pk-i, Pk-2, • •• After checking that the matrix elements (p\ph\p') and {h\pk\h') given by Eq. (|2.5| ) 
are consistent with the corresponding evolution equations (|2.4|) , one sees that the information 
contained in Eq. ( |2.4| ) concerns only the particle- hole and hole-particle matrix elements of pk- 

The corrections Ui (I = 1,2, ...) to the harmonic frequency u are determined by requiring 
that there is no resonant term with frequency ujq in the evolution equation of the particle-hole 
matrix elements of pi+i- We wish to emphasize that this procedure is sometimes inapplicable. 
Indeed, resonant terms with frequencies equal to 2co>o, 3uj , ... may appear in the construction 
of higher-order terms when these frequencies already belong to the continuous spectrum of the 
RPA matrix. However, one can circumvent this difficulty by generalizing the ansatz Q2.3p . The 
treatment of resonant terms is developped in Refs. ]16|, j24| . 

Once the periodic orbits are obtained, they can be used to construct energy spectra of N- 
body systems. The standard procedure is to perform a semiclassical quantization of these orbits 
p5| . This prescription consists in selecting those solutions for which the mean-field action / 
along a periodic orbit p(t), 

/ = £ [ T dt (Mmd t - %\Mt)) (2-6) 

is equal to an integer multiple of Plank's constant (or in atomic units, I = 27m). The summation 
runs over the quasi-periodic single-particle occupied states and Oh is the so-called Floquet- 
Lyapounov phase |16| . 

We now apply the method discussed above to the case of collective oscillations of metallic 
clusters (see Eqs. ( |1.2|) and (|L~6l ) ). In order to be able to discuss the two-phonon states, we 
perform the calculations up to the third order in the amplitude of the orbits. The explicit 
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expressions of the contributions Wk(k = 0, 1, 2, 3) to the mean-field Hamiltonian W appearing 
in Eq. fl2.4|) are 



Wo(r) = -~A + F (r) + / dV— L-p ( r ') + Vj(r) 
2 J \r — r'| 

Wi(r,t) = p^(r,t) + / dV- - Pl (r',t) 

dpo ^ | r — r 

(2.7) 

Ws(r,t) =^p 3 (r,t)+ /dV^-p 3 (r',t) + ^p 1 (r,t)p 2 (r,t) + i^p?(r,t) 
dpo J |r — r | dpo 6 dpg 

The quantity Fq(t) in this last equations is 

F (r) = -(-) 1/3 Po /3 W - 0.0333 Ln[l + 11.4(^p (r)) 1 / 3 ] (2.8) 

7T 3 

and the notation p k {r,t){k = 0,1,2,3) means (r|p fc (i)|r). 

The RPA matrix corresponding to the angular-momentum A can be written as (see, 

e.g., Ref. [0) : 



M w = I (2.9) 



where the matrix and are 



In equations fl2.10| ) we have introduced the notation a a = n a , l a where n a is the principal 
quantum number and l a the orbital angular momentum of the single particle static wave function 
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(r\(f) a ). The quantities e a and R aa {r) are respectively the single particle energies and the radial 
part of (r|0 a ). The coefficient / is given by 

l p lh A 




/£) =J(2l p + l)(2l h + l) 



while F (r) is already defined in Eq. (|2.8| ). The function Q appearing in these equations is : 

S%h<r) = I dr ' r ' 2 Ra P (r')R ah (r') 
where r< stands for Min(r, r') while r> = Max(r, r'). 

3 Results 

In this section we present some results obtained with the above formalism for monopole, dipole, 
quadrupole and octupole modes in charged Na^i- Calculations were performed on a lattice 
including 100 mesh points whose total size is 30 a.u., i.e., 15. 9A. In a first step we have 
solved the static equations. The values we have found for the binding energy, radius and single 
particle energies of Na^i are given in Table I. These values are in excellent agreement with 



those obtained by Catara, Chomaz and Van Giai |2T[ who use the same energy functional but 
a different numerical method. In a second step we have performed RPA calculations. The 
results for the energies of collective states of multipolarities L = 0, 1, 2, 3 and the corresponding 
percentages of the sum rule and energy weighted sum rule are collected in Table II. These 
values also agree well with those of reference p]J for collective modes. It can be noted from 
Table II that dipole and quadrupole modes exhibit a strong collective character while such 
is not the case for octupole and especially monopole modes. For the monopole mode one 
observes two neighbouring states with similar sum rule fractions. The monopole mode is a 
nearly unbound RPA state while other modes are bound. However, two-phonon states would 
belong to the continuum in the harmonic approximation. Transition densities of dipole and 
quadrupole states are displayed on Figs. 1 and 2 . The fact that these modes are bound is 
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reflected by the rapid decrease of these densities at large distance. For reference the static 
electron density is shown on the same scale (solid line). 

The observation of the dipole resonance in alkali-metal clusters including the charged Nd^i 
has been reported by several experimental groups, e.g., Brechignac et al |19[ and J. Pedersen 
et al [^(| (see also Ref. |26| for a review). A pronounced peak in the photoabsorption cross 
section is seen at ~ 2.7 eV. This observation agrees with the RPA calculations of Refs. |], [11 



These authors predict the dipole mode in Na^i at ~ 3.0 eV with « 80% of the sum rule and 
~ 83% of the energy weighted sum rule. These results also agree well with ours. 

In a third step we have determined the location of the not yet observed two-phonon states 
by performing a semiclassical quantization of our second order periodic orbits as prescribed by 
Eq. ( p.6|) . This prescription has been supplemented by an approximate angular momentum 



projection which generates a splitting of two-phonon states with different J's [ I6| , [24 |. For 
most states this procedure produces results which are stable with respect to the lattice size. 
The corresponding results are shown in Table III (a). One observes that small anharmonicities 
are found except for the J = two-octupole phonon state for which a 15 percent deviation is 
obtained. 

For the other states ( double monopole, J = 2 double dipole and J = 0,2,4 double 
quadrupole phonons) we were not able to calculate the energy shifts by the method described in 
Sec. II because resonant terms appear during the construction of second order terms. For exam- 
ple, the J = two quadrupole phonon state couples with a resonant mode of the RPA(L = 0) 
matrix so that the inhomogeneous linear equation Q2.4|) is singular in second order. To study 
these states we use the modified perturbative construction procedure presented in Ref. |L6| . As 
first order periodic orbit we take a superposition of the L-mode with energy Twol and the reso- 
nant RPA mode with energy huji r = 2%ujl instead of the ansatz (|2.3|) . With this modification, 
the semiclassical quantization procedure generates not only the splitting of the two-phonon 
states (n = 2) but also a splitting of each state n, J [|L6|1 . The results are shown in Table 111(b). 
We have checked that the results in this table change by less than 5 percent when the size of 
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the lattice is increased to 128 mesh points or decreased to 75 mesh points. One remarks first 
that the energy shifts are small and compatible with zero so that the anharmonicities of these 
states are negligible, second that for some states the energy shift acquires an imaginary part. 
We interpret this result as the signature of an instability of the corresponding orbit. 

A comparison with the results of Catara-Chomaz-Van Giai [^T[] for the modes in Table III is 
not straightforward. Indeed these authors find large anharmonicities arising from couplings be- 
tween states such as |(0 + , + )0 + > and |(2 + , 2 + )0 + > which are out of reach of our semiclassical 
quantization approach. Further studies of this question are now under way. 



4 Conclusion 

In this paper we have investigated the static and dynamic properties of the alkali-metal cluster 
iVaJi in the jellium model. Our results for the static and RPA calculations agree with those 
of reference |21| . We have determined the splitting of the two-phonon states of multipolarities 



L = 0, 1, 2, 3 in charged Na^i- This was achieved by using a perturbative method of constructing 
periodic orbits of the mean-field equations recently developped for giant resonances in nuclei 



16f . In most cases, we find small anharmonicities. To study the two-phonon states of the 
monopole and quadrupole modes we had to use a more general method which takes into account 
the coupling between the two-phonon state and the resonant modes fl6| . It turns out that these 
resonant terms generate anharmonicities which are so small as to be compatible with zero. For 
the J = 2 two-dipole, J = 2 two-quadrupole and J = 4 two-quadrupole states, we find an 
imaginary value for the energy shift with respect to the harmonic two-phonon energy. These 
values are however small and compatible with zero but are nevertheless a signature of an 
instability of the corresponding orbit. Such instabilities cannot be detected by other methods 
based on quasi boson expansions [JHJ]. 

Finally, we wish to emphasize that periodic orbits provide a powerful method to investigate 
collective energy spectra. This method has the advantage of being rather easy to implement. 
It was first found useful in the context of giant collective modes in nuclear physics [15|, |16| and 
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also in the context of baryonic physics to describe the Roper resonance |T7| . In both cases small 
anharmonic terms were found. This seems to be true also for giant collective modes in metallic 
agregates. This is a fortunate situation since expansions in the elongation of periodic orbits 
often break down, even for simple systems, because of the appearence of branches or bifurcations 
at small values of the elongation fll8 . Here in contrast the region of interest corresponding to 
one and two-phonon excitations lies well within the radius of convergence of the series. 
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Table captions 

TABLE I. Static properties of charged Na^i- The binding energy E LDA is in eV, the 
root-mean-square radius r in A and the single particle energies are in eV. 
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TABLE II. Low-lying RPA states of charged Nd^i- For each multipolarity L w we report 
the collective energy Huq (in eV) and the corresponding percentages of the sum rule (SR) and 
energy weighted sum rule (EWSR). 

L w hujo % SR % EWSR 
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TABLE III. Two-phonon states in charged Na^i- For each state are reported its spin and 
parity and the shift SE (in eV) with respect to the unperturbed two-phonon energy Eq = 2hu 
(in eV) (calculated with semiclassical quantization [See Eq. (|2.6|)1). 
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Figure captions 



FIG. 1. Static density po( r ) (^~ 3 ) (full line) and transition density pi(r) (A~ 3 ) (dashed 
line) corresponding to dipole state 1~ in the case of charged Na^i- 

FIG. 2. Static density po( r ) (^~ 3 ) (full line) and transition density pi(r) (A -3 ) (dashed 
line) corresponding to quadrupole state 2 + in the case of charged Na^i- 



